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Present knowledge of QCD n-point functions of Wilson lines at high energies is rather limited. In 
practical applications, it is therefore customary to factorize higher n-point functions into products 
of two-point functions (dipoles) which satisfy the BK evolution equation. We employ the JIMWLK 
formalism to derive explicit evolution equations for the 4- and 6-point functions of fundamental 
Wilson lines and show that if the Gaussian approximation is carried out before the rapidity evolution 
step is taken, then many leading order N c contributions are missed. Our evolution equations could 
specifically be used to improve calculations of forward dijet angular correlations, recently measured 
by the STAR collaboration in deuteron-gold collisions at the RHIC collider. Forward dijets in 
proton-proton collisions at the LHC probe QCD evolution at even smaller light-cone momentum 
fractions. Such correlations may provide insight into genuine differences between the JIMWLK and 
BK approaches. 

I. INTRODUCTION 

The Jalilian-Marian— Iancu— McLerran— Weigert— Leonidov— Kovner (JIMWLK) functional evolution equation de- 
scribes the energy dependence of n-point functions of Wilson lines at small light-cone momentum fraction x [TJ E] • 
These n-point functions appear in multi-particle production cross sections in hadronic (or heavy-ion) collisions. To 
date, our knowledge of the behavior of such n-point functions in QCD is very limited. Therefore, it is common to 
employ a Gaussian (and large- N c ) approximation which reduces these functions to powers of the two-point function, 
which at high transverse momentum corresponds to the well-known Balitsky-Fadin-Kuraev-Lipatov (BFKL) uninte- 
grated gluon distribution [3J. The evolution of the two-point function in the dipole approximation [1] is determined 
by an ordinary integro-differential equation known as the Balitsky-Kovchegov (BK) equation [TJ [5] which has been 
the subject of intense theoretical investigation in the past few years. 

Evolution with energy or rapidity y ~ log 1 jx occurs by (real or virtual) radiation of an additional gluon from a 
given n-point operator. Cross sections are related to expectation values of traces of such operators which project 
onto their physical matrix elements. If the expectation value of the n-point operator is split into dipoles before the 
radiation of the additional gluon (in order to perform the evolution step by means of the BK equation) then only the 
dipole from which the radiation emerged is allowed to split into two dipoles. On the other hand, if the evolution step 
is given by the JIMWLK equation then additional contributions arise, even at leading order in N c . This is illustrated 
in fig. 



FIG. 1: Left: BK-evolution of a set of dipoles: in an evolution step, emission of a gluon can only split the parent dipole. Right: 
JIMWLK evolution of the n-point function allows for additional contributions. 



To exhibit the differences between Gaussian+BK and JIMWLK evolution more clearly we focus on two specific 
n-point functions which appear in (forward) quark+gluon dijet production in pA collisions. Hadron-nucleus collisions 
are well suited for investigations of high gluon density QCD due to the fact that one can derive analytic relations 
for particle production when only the target is dense, while the projectile is dilute. Furthermore, hadron-nucleus 
collisions are free of the hot final state medium ( "quark-gluon plasma") produced in AA collisions which may affect 
the observables. Very recently, the STAR collaboration has presented two-hadron angular correlations in the forward 
rapidity region of deuteron-gold collisions which show a weakening/disappearance of the away side peak [6], in agree- 
ment with expectations from gluon saturation dynamics [T14L0] . The purpose of the present paper is to derive explicit 
and complete (in terms of iV c counting) evolution equations for the relevant n-point functions so that quantitative 
theoretical expectations could eventually be obtained. 

The cross section for production of a valence quark plus a gluon at forward rapidity in pA collisions was calculated 
in [11] (in momentum space) and in [12] (in coordinate space); we also refer the reader to refs. [13] . For completeness, 
we reproduce the expressions in the appendix. They involve expectation values of products of (up to) six Wilson lines 
in the fundamental representation (single- inclusive production involves only the 2— point function |14j). Explicitly, 
the following operators appear in the two parton (quark+gluon) production cross section 



4 (r, f : s) = tr V} t a V f t b [U s ] ab 



tr V s tr V f V% - — tr V? V f 



(1) 



and 



O a (r, f:s,s)= tr V r V} t a t b [U s Ul] ba 



tr V r V} V s V? tr V s V-] - — tr V r V} 



(2) 



where V (U) is a Wilson line in the fundamental (adjoint) representation and r, s etc. denote two-dimensional 
coordinates in the transverse plane. Here, we have used the identity 



jjab t b ^ V U a V 



(3) 



to relate matrices in the two representations. 

Due to the fact that explicit evolution equations for the above combinations of Wilson lines have not been derived 
so far, and therefore solutions to these equations are unknown, it is common to resort to the large N c and Gaussian 
approximation [10[ I12[ I15[ 116] . In this approximation, these two expectation values can be written as 

(0 4 (r,r:s)) ~ (0 2 (r - s)) (0 2 (s - r)> 
(0 6 (r,f :s,s)) ~ (0 2 (r-s)) (0 2 (f-s)) (0 2 (s-s)) + (0 2 (r-f)) (0 2 (s-s)) (0 2 (s-s)) (4) 

where (0 2 (r, f)) = (tr V r Vf). To keep track of factors of N c it is convenient to normalize expectation values of traces 
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of Wilson lines, i.e. the S matrices, as follows: 



S(r,f) ee Jj-<O a ) 



S 4 (r,f:s) ee -^{Oa) 

S 6 (r,f:s,s) ee ^^(Oe) . (5) 

For phenomenological applications of the Color Glass Condensate formalism to two-particle production, it is common 
to employ the approximations Q to factorize the higher point functions of Wilson lines into products of two point 
functions. This greatly simplifies applications as it allows one to write £4 and Sq in terms of the well known BK 
two point function. Then BK evolution of the two point function S is assumed to account for small-x evolution of 
the higher point functions S4 and Sq [TO]. We shall show below that indeed this procedure retains the leading- N c 
contribution to the 4-point function S4. On the other hand, we also show that it misses many leading- iV c contributions 
to the evolution of the 6-point function 5*6. 



II. EVOLUTION EQUATIONS FOR HIGHER POINT FUNCTIONS 

In this section we derive explicit evolution equations for the expectation values of O4 and Oq which appear in the 
two-particle production cross section. We then apply the large- N c and Gaussian approximations to our evolution 
equations and compare to the "naive" result Q where BK evolution for the two point functions is used on the right 
hand side of Q. 

We start from the JIMWLK evolution equation which determines the small- x evolution of any n-point function O 
from 

where 



„m = 1 f d " z (*-*)-(y-*) h+rfu _utu -utu] bd (7) 

We first focus on O4 defined in eq. . The second term (the two point function) from that equation evolves according 
to the JIMWLK equation 

> *> (r w (t. V}V, - ±* VJ V, *V,V!) . (8) 

This reduces to the BK equation in the Gaussian and large- N c approximations. On the other hand, the first term in 
O4 involving four fundamental Wilson lines satisfies 



d (tr V? V s tr V f y s t) = \(jd 2 x d 2 y ^ v % ^ tr V} V s tr V f 



dy 



(9) 



Using the explicit form of rftcy gi yen by (uk, we note that the derivative of r\ (w.r.t. a b x ) vanishes due to the color 
structure, except for the last term ~ U\\J y where one does need to differentiate U y . The calculation is straightforward 
but lengthy and involves repeated use of eq. ^ as well as the Fierz identity for the product of fundamental matrices 



[t a ]i] [t a }ki - 



(10) 
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Adding the evolution equation for the two point function, eq. (18]), we obtain 



(r - s) 2 



{f-s) 2 



1 



(r - s) 2 



(r — z) 2 (s — z) 2 (r — z) 2 {s — z) 2 
(f — s) 2 (r — f) 2 



tr Vl f V, tr V? Vf 



(r — z) 2 (s — z) 2 (r — z) 2 {s — z) 2 (r — z) 2 (r — z) 2 

(f-s) 2 



(r 2( S)2 _ , 2 tr Vj V z tr V? Vfc tr V, 



(r — z) 2 (s — z) 
(r - s) 2 

(r — z) 2 (s — z) 2 (f — z) 2 (s — z) 2 



(r — z) 2 (s — z)* 



tr V? K V? V? Vi + tr 7 r t y s v_ y; 



tr V; 1 K s tr V}_ V f tr V? V z 



(f-s) 5 



tr v; 1 y f 



1 ( _ ^T\2 

> L tr V< V tr \n V- 

N 2 (r-z) 2 (f-z) 2 K Vz Vz Vr 



(11) 



Note that the first and fourth line on the rhs could be combined into O4 by virtue of the identity ([!]); thus, all terms 
on the rhs feature at least as many Wilson lines as the original operator. 

Equation ( 11 ) describes the small- a; (rapidity) evolution of the 4-point function O4. It has the intuitive interpretation 
of the scattering matrix for two quark-anti-quark dipoles evolving with rapidity while (multiply) scattering from the 
target. It is worth noting that this equation remains finite when the internal integration variable z approaches any of 
the external coordinates r, f, s. 

We now apply the large N c and Gaussian approximations to eq. (11). The leading- N c contributions on the rhs 
originate from the first and third lines of ( 11 ). Written in terms of the scattering matrix S4, 

(r-s) 2 



dy 



Sa 



N c a s 
' 2tt 2 



d l z { S{s - f) 



+ S(r-s) 



Using the BK equation for the 2-point function S, 



dy 



S(r-s) 



N c a s 
2tt 2 



d 2 z- 



(r — z) 2 (s — z) 2 

{f-s) 2 
(f — z) 2 (s — z) 2 



(r - s) 2 



(r — z) 2 (s — z) 2 



S(r - s) - S(r - z) S(z - s) 
S(f — s) — S(f — z) S(z — s) 

S(r — s) — S(r — z) S(z — s) 



we can rewrite eq. ( 12 ) as 



— SAr, r : s) ~ Sis — r) — S(r — s) + Sir — s) -^-S(s — r) 
dy dy dy 



(12) 



(13) 



(14) 



Thus, it is evident that the factorized form Q for the 4-point function (O4) correctly captures the leading- N c 
contributions. In passing, we note that to arrive at eq. (14) we have dropped 12 terms of order 1/N 2 . 

Next, we derive the evolution equation for the 6-point function Og from the JIMWLK equation The procedure 
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is straightforward but very tedious, here we just quote the final result: 

d 2 z 



(r-sf 



(r — rY 



(f-s) 2 



+ 



(r — z) 2 (s — z) 2 (r — z) 2 (r — z) 2 (r — z) 2 (s — z) 2 (s — z) 2 (s — z) 2 
(r — s) 2 (r — s) 2 (s — s) 2 



trVrVlVsVj trV s Vl 



-[- 



+2 

+ 



(r — z) 2 (s — z) 2 (r — z) 2 (s — z) 2 (s — z) 2 (s — z) 2 
(f-s) 2 (f-s) 2 (s-s) 2 



+ 



r — z) 2 (s — z) 2 (r — z) 2 (s — z) 2 (s — z) 2 (s — z) 2 ■ 
(f-s) 2 (f-s) 2 (r - s) 2 



r — z) 2 (s — z) 2 (r — z) 2 (s — z) 2 (r — z) 2 (s — z) 2 (r — z) 2 (s — z) 2 

(r - s) 2 



]trv r viv s v2v t viv a v2 
trv r viv x viv.v2v s v2 

}trV r V}V a vlv z v}V- s Vl 



(r - sf 



(f-s) 2 



... i \ i > \ ■ - „ < \ i , i \ ■-> .... j - ! \ '"■ - ' 

(r - f) 2 



(r - r) 2 



+ 



(r - s) 2 



- z) 2 (s — z) 2 
(f - s) 2 
~z) 2 (s-z) 2 

(8-8)* 



+ 



s — z) 2 (s — z) 2 
(f-s) 2 





(f - s) 2 




(r - s) 2 




(f 


-z) 2 (s- 


z) 2 


(r — z) 2 (s — 


z) 2 




(r-s) 2 




(s - f) 2 




(r 


-z) 2 (s- 


z) 2 


(s — z) 2 (f — 


z) 2 




(f-s) 2 




(r - s) 2 




(f 


-z) 2 (s- 


z) 2 


(r — z) 2 (s — 


z) 2 




(s-s) 2 




(r - s) 2 




(s 


-z) 2 (s- 


z) 2 


(r — z) 2 (s — 


z) 2 




(S-S)> 




(f-s) 2 




(s 


-z) 2 (s- 


z) 2 


(f — z) 2 (s — 


z) 2 


V r VlV s v2 


tr V s Vt tr V z Vi 






(f - s) 2 




(r - s) 2 




(f 


-z) 2 (s- 


z) 2 


(r — z) 2 (s — 


z) 2 




(r — i 


■? 


(r — 7 





tr V r Vl V s Vf V s V} 



}trV z V}V- s V} trVrV? tvV s V^ 

] tr V r Vl V, V* tr V z V} tr V s V^ s 

]\xV T V}V s V} tvV z V} tvV a vl 

]\xV T V}V x V} trVsV} trVsVl 



+ 



(r - s) 2 



tr V r K f tr V} V s tr V s V 



^- S ) 2 ( r ~ 5 ") 2 + ( r - f ) 2 + ( S -^ 2 l trF y_i tr y_ y t toVVl 

(f-z) 2 (s-z) 2 (r-z) 2 (s-z) 2 ^ (r-z) 2 (f-z) 2 ^ (s - z) 2 (s - z) 2 1 r r s s 



(r - s) 2 (r - s) 2 



(f-s) 2 



(r - s) 2 



(r — r) 



r — z) 2 (s — z) 2 (r — z) 2 (s — z) 2 (r — z) 2 (s — z) 2 (r — z) 2 (s — z) 2 (r — z) 2 (r — z) 2 

(8-8) 2 



(r — r) 2 



N 2 (r-z) 2 (f-z) 2 



trV r V z trV^yt 



(15) 



We have checked that this equation remains finite when the internal integration variable z approaches any of the 
external coordinates r, f , s, s. Note the appearance of lower point functions on the right hand side of the equation. 
Unlike the similar terms in the evolution of O4 these can not be combined into the original operator Og. 

We may now employ the Gaussian approximation on the rhs of the equation above to exhibit the lcading-iV c 
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contributions (from lines 1 and 7 — 13). The equation reduces to 



dy 



Sri 



N c a s 
" (2tt) 2 



(r — r) 



d 2 z 



(r-s) 2 



(r — r)^ 



(f - s) : 



(s - s) : 



(r — z) 2 (s — z) 2 (r — z) 2 (r — z) 2 (r — z) 2 (s — z) 2 (s — z) 2 (s — z) 2 
S(r - f) 5(a - s) 5(« - s) + S(r - s) S(f - s) S(s - s) 



(r - s) 2 



(s - ff 



(r — z) 2 (r — z) 2 (r — z) 2 (s — z) 2 (s — z) 2 (r — z) 2 
(r — s) 2 (r — f) 2 (f — s) 2 



S{z - r) S(s - s) + S{z - s) S(f - s) S(r - z) S{s - s) 



(r 


-z) 2 (s- 


z) 2 




(r - s) 2 




(r 


— z) 2 (s — 


z) 2 




(f - s) 2 




(f 


-z) 2 (s- 


z) 2 




( S - S -) 2 




(s 


-z) 2 ( S "- 


z) 2 




(f - s) 2 




(f 


— z) 2 (,s — 


z) 2 



( S - S -) 2 



(r - s) 2 



( S - S -) 2 



(f- s) 2 



z)5(s- 


5) - 


fS(r- 


s)S'(z 


f)S(s- 




h5(r- 


«)5(f 


f) S(z - 


»)- 


h5(r- 


s)5(f 



5(r - r) S(s - s) + S(r - s) S(f - s) J S(s - z) S{z - s) 
(f — s) 2 (r — s) 2 (r — s) 2 



r — z) 2 (s — z) 2 (r — z) 2 (s — z) 2 (r — z) 2 (s — z) 2 

(f— s) 2 (r — s) 2 (r — r) 2 (s — s) 2 

(f — z) 2 {s — z) 2 (r — z) 2 (s — z) 2 (r — z) 2 (f— z) 2 (s — z) 2 (,s — z) 2 



5(r — s) — s) S(s — s) 
S(r — f) 5(s — s) 5(s — s) 



There are many more terms in this equation than obtained by differentiating eq. Q. To make this more clear, we 
write eq. ( 16 ) in the form 



d a i - -\ d 
dy dy 

N c a 8 f 2 

— S(s-s)J dz 

(r-s) 2 ( 
(r — z) 2 (s — z) 2 



5(r - s) S(S - f) S{s - s) + S{r - f) S(s - s) S(s - s) 



+ 



+ 



(r — r) 2 
(r — z) 2 (f — z) 2 



(r — s) 2 
(r — z) 2 (s — z) 2 

(r-s) 2 
(r — z) 2 (s — z) 2 



s 



s 



( 

(f - s) 2 


(r 


— z) 2 (s — z) 2 


?(r 


- f) - 5(z - 




(s-s) 2 


(s 


— z) 2 (s — z) 2 


>(r 


- s) - S(z - 




(r — f) 2 


(r 


— z) 2 (f — z) 2 




(r - s) 2 


(r 


— z) 2 (s — z) 2 



(f - ,s) s 



S(s - s) 



(f- s) 2 



S(f - s) 



( S -s) 2 



S ") 2 



(r - s) 5 



(f - s) 2 



(f - s) 2 



(f - sf 



S(r — s) S(z — s) S(z — r) 
S(r — f) S(s — z) S(z — s) 



(17) 



The terms in the bracket (first line) could also be obtained from the Gaussian approximation Q to the 6-point 
function combined with BK evolution of the dipoles 1 . All other terms would be missed. The structure of the extra 
terms suggests why a naive Gaussian approximation fails; to see this, consider Oq from eq. ([2]). The first term is a 
product of two traces, a trace of 4 Wilson lines times a trace of two Wilson lines. Emission of a gluon between a different 



(16) 



1 The second term in the bracket in the first line seems to have been omitted in ref. 1 1 21 . 
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quark-anti-quark pair is missed by the Gaussian approximation. In other words, if the Gaussian approximation is 
applied to the trace of four Wilson lines in eq. ^ before the rapidity evolution step, then due to color neutrality, 
a radiated gluon can not end up in another dipole. On the other hand, if one performs the rapidity evolution step 
of the full 4-point function then the emitted gluon can end up anywhere between other quark and anti-quark pairs, 
including those which are not the "parents" of the radiated gluon. These are precisely the contributions which are 
missed by the Gaussian approximation in the case of the 6-point function 2 . It is easy to see that this is not possible for 
the 4-point function ([T]) since the traces involve only two Wilson lines, and only one Gaussian contraction is possible 
(at leading order in Ay. This also means that our findings here do not affect fully inclusive observables such as DIS 
structure functions F2 and Fl, or single inclusive particle production in pA collisions, since those involve the two 
point function; its evolution is only sensitive to the 4-point function which does not receive any leading N c corrections 
from JIMWLK (as compared to BK). 

In summary, we have derived explicit evolution equations for the n-point functions that appear in forward dijet 
angular correlations in pA collisions. We find that factorizing these n-point functions into dipoles before performing 
evolution in rapidity misses many leading- A" c contributions; higher multipole operators obey different evolution equa- 
tions which can not be reduced to BK evolution of dipoles [2T]. Moreover, a rather large number of A^-suppressed 
terms arises in the full JIMWLK evolution equations for higher n-point functions which may give substantial numerical 
contributions, especially when n > N c . 

The results presented here underscore the importance of rigorous solutions of the small- a; evolution of higher point 
functions of Wilson lines. These could be obtained numerically via lattice-gauge theory techniques along the lines of 
ref. [EH HO] where the small x evolution of the two point function has been studied. For the two-point function, those 
authors found only very minor differences between the JIMWLK and BK (Gaussian + leading N c approximation) 
evolution equations. One may expect that due to the many leading- A^ terms missed by the Gaussian approximation 
the differences between the JIMWLK and BK evolution of the 6-point function should be substantial (see also [TBI |22J 
for other observables where possible differences between JIMWLK and BK were investigated). Once (numerical) 
solutions to these evolution equations become available, they could also be used to improve present calculations [7rtT0] 
for dijet production and angular correlations. 



Appendix A: Forward q + g production in pA collisions 

Here, we reproduce the expression for forward q + g production (with comparable rapidities) in a valence-quark 
nucleus collision from ref. |12| . This expression is to be convoluted with valence quark distribution functions of a 
proton (or deuteron) and with parton — > hadron fragmentation functions in order to obtain a physical cross section, 
see refs. [101112]. ' 

The qA — > qgX hard scattering cross section is given by 

d(J qA^ qg X ^ na CjpS(t) + -k + -a + ) f ^ r ik±-(x'-x)+i( q± -p ± )-(b'-b) 

a?ka?q ~ a S^F0[p ft q )J ^ ( ^ )2 ^ { ^ )2 
J2€p(P,k+,x' -b')^(p,k+,x-b) 

Aa/3 

[S 6 (b, x, b', x) - S 4 (b, x, b' + z{x' - b')) - S 4 {b + z{x - b),x', b') + S(b + z(x -b),b' + z{x' - b'))] . (Al) 

Here, z — k + /p + and (f>*cj(p,k + ,x) denotes the amplitude of the \qg)o Fock state component in the wave function of a 
dressed quark to leading order in a s ; the explicit expression is given in ref. [T2]. Lastly, the various n-point functions 
(target averages) are given by 

S 6 (b, x, b', x') = ^V(tr V b V b U d t c [U x Ul,] cd ) , (A2) 
5 4 (6, x, b 1 ) = (tr Vj,t c V b t d Uf) , (A3) 

S(b,b') = ~{tvV b V b \) . (A4) 



In this respect, the effects discussed here go beyond the correlations seen in dipole chains of the form tr V x Yj tr V r Vf tr V s Vj ■ ■ ■ [T7] . 
For such operators there is no cross-dipole emission as in fig. [T] 
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